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We study both analytically and numerically the complexity of the adiabatic quantum evolution 
algorithm applied to random instances of combinatorial optimization problems. We use as an exam- 
ple the NP-complete set partition problem and obtain an asymptotic expression for the minimal gap 
separating the ground and exited states of a system during the execution of the algorithm. We show 
that for computationally hard problem instances the size of the minimal gap scales exponentially 
with the problem size. This result is in qualitative agreement with the direct numerical simulation of 
the algorithm for small instances of the set partition problem. We describe the statistical properties 
of the optimization problem that are responsible for the exponential behaviour of the algorithm. 

PACS numbers: Valid PACS appear here 


Since the discovery by Shor [1] nearly a decade ago 
of a quantum algorithm for efficient integer factorization 
there has been a rapidly growing interest in the devel- 
opment of new quantum algorithms capable of solving 
computational problems that are practically intractable 
on classical computers. Perhaps the most notable exam- 
ple of such problems is that of combinatorial optimization 
(COP). In the simplest case the task in COP is to mini- 
mize the cost function (” energy” ) E z defined on a set of 
2" binary strings z = {zi, . . . , z n j Zj = 0, 1, each contain- 
ing n bits. In quantum computation this cost function 
corresponds to a Hamiltonian Hp 

h p = j 2 e * |z><z| (i) 

z 

|z) = |zi)i <g> \zi)l ® ® | Z n ) n - 

where the indices Zj =0,1 and the summation is over 2" 
states |z) forming the computational basis of a quantum 
computer with n qubits. State \zj)j of the j-th qubit 
is an eigenstate of the Pauli matrix a z with eigenvalue 
Sj — 1 — 2 Zj ( Sj = ±1). It is clear from the above that 
the ground state of Hp encodes the solution to the COP 
with cost function E z . 

COPs have a direct analogy in physics, related to find- 
ing ground states of classical spin glass models. In the 
example above bits Zj correspond to Ising spins Sj. The 
connection between the properties of frustrated disor- 
dered systems and the structure of the solution space 
of complex COPs has been noted first by Fu and An- 
derson [2]. It has been recognized [3] that many of the 
spin glass models are in almost one-to-one correspon- 
dence with a number of COPs from theoretical computer 
science that form a so-called NP-complete class [4]. This 
class contains hundreds of the most common computa- 
tionally hard problems encountered in practice, such as 
constraint satisfaction, traveling salesmen, integer pro- 
gramming, and others. NP-complete problems are char- 
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acterized in the worst cases by exponential scaling of the 
running time or memory requirements with the prob- 
lem size n. A special property of the class is that any 
NP-complete problem can be converted into any other 
NP-complete problem in polynomial time on a classical 
computer; therefore, it is sufficient to find a deterministic 
algorithm that can be guaranteed to solve all instances 
of just one of the NP-complete problems within a poly- 
nomial time bound. However it is widely believed that 
such an algorithm does not exist on a classical computer. 
Weather it exists on a quantum computer is one of the 
central open questions. Ultimately one can expect that 
the behavior of the new quantum algorithms for COPs 
and their complexity will be closely related to the prop- 
erties of quantum spin glasses. 

Farhi and co-workers suggested recently a new quan- 
tum algorithm for solving combinatorial optimization 
problems which is based on the properties of quantum 
adiabatic evolution [5]. Running of the algorithm for 
several NP-complete problems has been simulated on a 
classical computer using a large number of randomly gen- 
erated problem instances that are believed to be com- 
putationally hard for classical algorithms [6-8]. Results 
of these numerical simulations for relatively small size 
of the problem instances ( n < 20) suggest a quadratic 
scaling law of the run time of the quantum adiabatic 
algorithm with n. In [5, 9] special symmetric cases of 
COP were considered where symmetry of the problem al- 
lowed to describe the true asymptotic behavior (n -» oo) 
of the algorithm. In certain examples considered in [9] 
quantum adiabatic algorithm finds the solution in time 
polynomial in n while simulated annealing requires ex- 
ponential time. However, so far there are no analytical 
results on the asymptotic behavior of the algorithm for 
the general case of randomly generated hard instances of 
NP-complete problems. 

In what follows we derive the asymptotic complexity 
of the quantum adiabatic algorithm for the Set Par- 
tition Problem (SPP). It is one of the six basic NP- 
complete problems that are at the heart of the theory 
of NP-completeness [4]. It can be formulated as com- 
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binatorial optimization problem: Given a sequence 

of positive numbers {do, «i ..... . a n ) find a partition, i.e. 

two disjoint subsets A\ and A 2 such that the residue 


E = 



is minimized. Due to the ob- 


vious symmetry of the problem it is always possible 
to fix the assignment of one of the numbers, say do, 
to a subset A\. In SPP we search for the bit strings 
z = {zi, . . . ,z n j (or corresponding Ising spin configura- 
tions S = {Si, . . . , S n j) that minimize the energy or cost 
function E z 


n 

Ez. = |0s| , Os = do + Sj = 1 — 2 Zj, (2) 

j 1 

Here Sj = 1 ( Zj = 0) if aj € A and Sj = — 1 
(zj = 1) if dj € A'. The minimum partition can also 
be viewed as a ground state of the Ising spin glass, — 0|, 
with the Mattis-like (global) antiferromagnetic coupling, 
Jij — d / dj . 

SPP also has many practical applications including 
multiprocessor scheduling [10], cryptography [11], and 
others. The interest in SPP also stems from the remark- 
able failure of simulated annealing to find good solutions, 
as compared with the solutions found by deterministic 
heuristics [12]. The apparent reason for this failure is 
due to the existence of order 2" local minima whose en- 
ergies are of the order of 1 /n [13] which undermines the 
usual strategy of exploring the space of configurations S 
through single spin flips. 

The computational complexity of random instances of 
SPP depends on the number of bits b needed to encode 
the numbers dj. Numerical simulations show [14, 15] that 
for independent, identically distributed (i.i.d.) random 
b-bit numbers dj, the solution time grows exponentially 
with n for n < b and polynomially for n > b. The tran- 
sition from the ” hard” to computationally ” easy” phases 
has features somewhat similar to phase transitions in 
physical systems [16]. Transitions of this kind were ob- 
served in various NP-complete problems [17] and can be 
analyzed in the framework of statistical mechanics. In 
what follows we concentrate on the computationally hard 
regime n <C b. 

We now consider the distribution of the partition ener- 
gies E z . It was derived in [18] using statistical averaging 
over the assemble of instances of SPP with i.i.d. num- 
bers dj. Bearing in mind our eventual goal of deriving the 
complexity of the quantum algorithm for a given random 
instance of SPP we do not perform such an averaging. 
Instead we work with a given set of {dj} and introduce 
a coarse-grained “density of states” 

P(E) = - E a )). (3) 

Z 

Here the sum is over 2" bit-strings z; partition energy 
E z is defined in (2) and (• • • ) denotes averaging over the 
interval of the partition values (E—r]/2,E+r]/2) with the 
window size chosen self-consistently, r] l/p(E). Using 


(2) we can rewrite (3) in the form 


p(E) = J q ds f (-y) cos(Es), (4) 

where I(s) = {IjLo cos ( a j s ) an ^ f( x ) = sin(x)/x is a 
window function that imposes a cut-off in the integral 
(4) at s ~ 2/r] <C p(E). For large n this integral can 
be evaluated using the steepest descent method. To find 
the saddle points we shall assume that the &-bit num- 
bers dj are distributed inside of the interval (0, 1) and 
are the integer multiples of 2~ b , the smallest number 
that can be represented with available number of bits. 
We note that for large n the function I(s) has sharp 
maxima (minima) with width ~ n -1 / 2 at the points 
St = kn2 b , k = 0,1,...; |I(s/j,)| = 1. The sum over 
these saddle-points was evaluated by Mertens [16] who 
analyzed the partition function in SPP. In our case only 
one saddle point at s = 0 contributes to the integral (4). 
Indeed, it will be seen below that the window size 2/r] 
can be chosen to obey the conditions 1 <C n 1 / 2 /r] <C 2". 
Therefore in the high-precision case ( b n) saddle- 
points Sfc with k > 0 lie far outside the window and 
their contributions can be neglected [? ]. On the other 
hand the window function f(x) can be replaced by unity 
while computing the contribution from the saddle-point 
at s = 0. Finally we obtain for E <C n 


2 n+1 0(E) 
\J2 7T ( 


E 2 


P(E)= ==~ 6xp 1 2,7^ „ 

<r 2 — 1 “J + Ofn^ 1 ), E -C n. 

U 3 = 1 


+ 0(n~^ 2 ) 


( 5 ) 


Here O denote a step function. 0(x)=l, for x > 0 and 
0(x)=O for x < 0. The coarse-grained state density in 
Eq.(5) depends on the set of {«j} through a single self- 
averaging quantity <r 2 (cf. [16]). For E that are not 
too large (E < \/n) the average separation between the 
individual partition energies E m i n ~ ■ s /n2~ n . This justi- 
fies the choice of the window r] above that corresponds 
to coarse- graining over many individual level separations. 
For a typical set of numbers {aj} there are only two scales 
present in the distribution of partition energies, one is a 
“microscopic” scale given by E m j n and another is given 
by the mean partition energy (E) = a(2n/n) 1 P 2 . 

Using the same approach one can derive the distri- 
bution p(O) of Ds- For |fis| “C n it is Gaussian, 
/5(±f2) = 1/2/9(101) and can be understood in terms of 
a random walk with coordinate 0 (cf. (2)). The walk 
begins at the origin, 0 = 0, and makes n steps. At the 
j-th step O moves to the right or to the left by “distance” 
2 Oj if Sj = 1 or Sj = — 1, respectively. In the asymp- 
totic limit of large n the result (5) corresponds to equal 
probabilities of right and left moves and the distribution 
of step lengths coinciding with that in the set of numbers 

{2 cij}. 

We now consider the “conditional” density of the parti- 
tion energies p r (E'\E) = p r (E',E)/p(E) where p r (E'\E) 
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is a “joint” density defined as follows 

p r (E',E) = Y l PrAE')S(E-E z ) (6) 

Z 

PrAE') = (J2 6 ( E ' ~ E *') A r, £>(■',■)) (7) 

z' 

Here is a Kronecker delta. Function D(z,z') above 
computes the number of bits that take different values 
in the bit-strings z and z', it is a so-called Hamming 
distance between the strings 

n 

D ( z , z ') = z 'i + C 1 “ ( 8 ) 

j=i 

Conditional density p r ,z(E') in (7) determines the dis- 
tribution of energies E z i for the strings z' obtained 
from z by flipping r bits. Similar to Eq. (3), {•••) in 
(6), (7) denotes averaging over the small energy interval 
(E — r)/2,E + rj/2) that, however, includes many indi- 
vidual partition levels. For r = 1 there are exactly n 
strings z' on a Hamming distance r =1 from the string 
z. Partition energies corresponding to these strings equal 
|fis — 2aj5j |, 1 < j < n (cf. (2)). One can show that for 
E z < 2 the distribution of these energies is nearly uni- 
form in the range (0, 2 — E z ) due to the fact that spins 
Sj are nearly equally distributed between ±1 values for 
configurations with |Hs| “C n. Similar arguments apply 
for the case r = n — 1. Using a coarse-graining window 
such that n~ l « g <C 1 one can obtain in the range of 
energies E',E Z < 1 the following result: 

p r (E'\E z )np rtZ (E')=n/2 + 0(l) 
r = 1 , n — 1 , n » 1 . ( 9 ) 


For large r,n — r 1 one can calculate the condi- 
tional densities (6) for a given set of {aj} by evalua- 
tion of the appropriate integrals with the steepest descent 
method, in a similar manner to the derivation of of p(E) 
in Eqs.(3)-(5). Here we will skip the derivation and only 
provide the result that we found to be in a very good 
agreement with the numerical simulations of SPP even 
for relatively small n ~ 20 — 30. State density functions 
in (6) are given each by a sum of Gaussians for E,E' > 0 
and have a broad maximum at the origin. Near the max- 
imum we have: 


Pr(E'\E z ) W Pr,z(E') W 

C(x) = (4x(l -x))“ 1/2 , 


C (r/n) fn\ 
a^/nn/2 v ) 
r,n-r>l. 


(10) 


This result applies in a range E 1 , E z <C [nc(r /n)] 1 ^ 2 . We 
note that conditional density (7) is uniform there and 
does not depend on E' and z after the coarse-graining 
over the small interval r) of energies E z < in (10). Window 
size ??>(") X . 

Eqs. (9), (10) suggest the structure of the low-lying par- 
tition configurations. Bit strings z' that are on a Ham- 
ming distance r from a given string z correspond to a sub- 
set of partition energies with average spacing ~ (") 


which does not depend on z in the range E Z ,E Z * < 1 
for all r. The greater the Hamming distance between 
the strings - the closer their partition energies can be to 
each other - and vise versa. This is a basic property of 
SPP that gives rise to an exponentially large number of 
local minima in the relatively narrow range of energies 
0(1). In a typical instance of random SPP correlations 
between the Hamming distances to a given string (e.g., 
the solution to SPP) and partition energies essentially 
vanishes. 

In the quantum adiabatic algorithm [5] one specifies 
the time-dependent Hamiltonian H (t) = H(t/T) 

H(t) = (1-t)V+tH p , (11) 

where r = t/T is dimensionless “time”. This Hamil- 
tonian guides the quantum evolution of the state 
vector | ip(t)) according to the Schrodinger equation 
id\ip(t))dt = H(t)\ip(t)) from t = 0 to t = T, the run 
time of the algorithm. Hp is the “problem” Hamiltonian 
given in (1). V is a “driver” Hamiltonian, that is de- 
signed to cause the transitions between the eigenstates 
of Hp. In this algorithm one prepares the initial state of 
the system ip(0) to be a ground state of H( 0) = V. In 
the simplest case 

n— 1 

V = -Y,°l IV’(O)) = 2“" /2 Xl z ); (12) 

j = 0 Z 

where al is a Pauli matrix for j-th qubit. Consider in- 
stantaneous eigenstates \<t> a (j)) of H(r) with energies 
A a (r) arranged in nondecreasing order at any value of 
re (0,1) 


H \4> n )=\ n \4>u), V = 0,1,..., 2" -1. (13) 

Provided the value of T is large enough and there is a 
finite gap for all t e (0, T) between the ground and ex- 
ited state energies, g(r) = Ai(r) — Ao(r) > 0, quantum 
evolution is adiabatic and the state of the system | ip(t)) 
stays close to an instantaneous ground state, \(j)o(t/T)) 
(up to a phase factor). Because H(T) = Hp the final 
state \ip(T)) is close to the ground state \4>o(t = 1)) of 
the problem Hamiltonian. Therefore a measurement per- 
formed on the quantum computer at t = T (r = 1) will 
find one of the solutions of COP nearly with certainty. 
In the Landau-Zener picture a quantum transition away 
from the adiabatic ground state occurs most likely in the 
vicinity of the point r* € (0,1) where the energy gap 
g{r) reaches its minimum (avoided-crossing region). The 
probability of the transition, 1 — |(^>(T)|0 o (t = 1))| 2 , is 
small provided that 

T>> A^’ 5 min = min [Ai(r) - A 0 (r)] (14) 

Here i I = \{4>i |-^|^o)|r=r* • Even in the cases that are 
more complex than the Landau-Zener “2- level” picture 
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of nonadiabatic transitions the analysis of minimum gap 
Smin can provide an important test for the complexity of 
the quantum adiabatic algorithm. E.g., if one can show 
that for a random COP ry min decreases exponentially fast 
with the problem size n while |i?oil does not decreases 
with n that would be a strong indication for exponential 
complexity of the algorithm in such problem. 

As suggested in [5] the quantum adiabatic algorithm 
can be recast within the conventional quantum comput- 
ing paradigm using the technique introduced by Lloyd 
[19]. Continuous-time quantum evolution can be approx- 
imated by a time-ordered product of unitary operators, 
e -t{i-T k )vs e -tT k H P 5 ^ corresponding to small time inter- 
vals {tk,tk +5). Operator e~ l (1-T 0 V 6 typically corre- 
sponds to a sequence of l(2)-qubit gates (cf. (12)). Op- 
erator e ~ lT k HpS i s diagonal in the computational basis |z) 
and corresponds to phase rotations by angles E Z A. Since 
in the case n <C b, the average separation between the 
neighboring values of E z is l/p(E) = 0(2~ n ), the quan- 
tum device would need to support a very high precision 
of physical parameters (like external fields, etc) to control 
small differences in phases, 0(2~ n ). Since this precision 
scales with n exponentially it would strongly restrict the 
size of an instance of SPP that could be solved on such 
quantum computer (this technical restriction is generic 
for COPs that involve high-precision numbers). To avoid 
this restriction we define a new cost function that takes 
values on a scale which is logarithmic in l/p(E). I.e., we 
chose the cost to be an oracle-type function £ z that takes 
a set of M + 1 integer values £& = — M + k, 0 < k < M 

£ z = -k + M, for u)k < E z < u>k+ i, (15) 

w* = ( 2 fc -l)w 0 , k = 0,...,M. (16) 

We set £ z = sm = 0 for all states with E z > % where 
com is a control parameter chosen sufficiently large so 
that M = 0(n). The space of 2" states |z) is divided 
into M + 1 subspaces £*,, each determined by Eqs. (15) 
and (16) for a given k 

M 

= E i z >< z i- ( 1? ) 

fc=o 

Note that subspace Co contains the solution(s) to the 
SPP. Dimension do of Co is controlled by the value of 
cuo in (16) which is another control parameter of the 
algorithm. We set u>o = K/p(0) where the integer 
K m do 1 is independent of n and determines how 
many times in average one needs to repeat the quantum 
algorithm in order to obtain the solution to SPP with 
probability close to 1. 

Operator Hp projects any state | ip) onto the states 
with partition energies in the range 0 < E z < com- We 
assume that com <C (E) so that the density function (5) is 
nearly uniform in this range. Therefore the dimensions of 
the subspaces Ck grow exponentially with k : dk = do 2 k 
for k < M. This simplification does not affect the com- 
plexity of a quantum algorithm that spends most of the 


time in “annealing” the system to much smaller partition 
energies, % > E z ~ E min = 0(n 1 ^ 2 2~ n ). 

We now consider a stationary Schrodinger equation 
(13) and obtain the minimum gap g m i n (14) in the asymp- 
totic limit n -» oo. To proceed we need to introduce a 
new basis of states |x) = |xi)i <g> |x 2 )2 ® ■■■ ® |x„)„ 
where state | Xj)j is an eigenstate of the Pauli matrix a x 
fo the j-th qubit with eigenvalue 1 — 2 Xj = ±1. Driver 
Hamiltonian V can be written in the following form 

n 

v =Y J v m i m , i m = l x >(4 ( 18 ) 

m= 0 x\-\ \-Xn=m 

For a particular case given in Eq. (12) we have V m = 
2m — n. Matrix elements of I m in a basis of states |z) 
depend only on the Hamming distance D( z, z') between 
the strings z and z' 


'TUI TUI 

■^z,z' *D( z,z') 


(19) 


i? = 2 -EE(": r )Q(- 1 ) p A ™.^ ( 2 °) 

We now rewrite Eq. (13) in the form 

W= X^aV Hp ^ )j a = a ( T ) =1 ~ T - ( 21 ) 

(we temporarily drop the subscript p indicating the num- 
ber of a quantum state and the argument r in and A) . 
^From (15)-(21) we obtain the equation for the ampli- 
tudes \t z = (ft|z) in terms of the coefficients I™ 


2~ n 

“ T 7 T l£d(z,z') £z' 3 2 G £k 


A — aVo 


z'^z 


= S*(A) = 1 -Tj 0 (X)e k , $ = (23) 


J r = J T (X) = Y, 


Ul= 1 


A — aVm 


0 < r < n. 


We now make a key observation that the (coarse-grained) 
distribution of matrix elements Jd(z,z') for a given z is 
determined by J r { A) and conditional density p r%z >{E) (6). 
We compute a cumulative quantity 


z^z 7 

J2 Mz,z')W » F kW + f.,kW (25) 

zee*, 

4(A) = I®, 8( A) = dr C,{r /n) (") J r (A) (26) 

where p = Dm/ yjnna 2 /2 and f z ,k( A) is a small correc- 
tion that will be described below. Function s( A) can be 
evaluated using the explicit form of J r (X) that decays 
rapidly with r. E.g., for |A — aVo| <C 1 one has 


Jr (A) = 


n 


— 1 n—r 2 —n f 71 ) 
\m+r/ 


E 

Ul= 1 


m 


(lnr + 7) . (27) 
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(7 is an Euler’s constant) and s( A) w — ln2/(2a). Ac- 
cording to (27), J r w (n/2 - r) _1 (") 1 for n/2 - r » 1, 
and the integrand in s(A) is a smooth function of r. This 
corresponds to the interesting fact that contributions to 
the sum in(25) from states z with D( z,z') = r <C n are 
small (0((r /n) 1 / 2 ). In fact, the dominant contribution 
is of “entropic” nature: it comes from the large number 
of states z € corresponding to Hamming distances, 
1 C r < n/2, and contributing each a small weight J T . 
The energies of these states are distributed, for a given r, 
according to p T}Z '(E) m p r . and this explains why F*( A) 
does not depend on z' (cf. (10)). 

We note that the entropic term Fk scales down with k 
exponentially fast (oc 2 k ~ M ), proportionally to the size of 
the window LOk+i — w k of partition energies corresponding 
to z £ Ck (cf. (16)). When this window becomes too 
narrow the contributions to the sum in (25) from the 
states z with r = D( z,z') ~ 1 are no longer described 
by the density of states p r , due to the discreteness of 
the partition energy spectrum E z . E.g., for a given z' in 
this case there could be either one or none of the states 
z € £k with D( z,z') = 1, depending on the choice of 
z'. Contributions to the sum in (25) from these states 
strongly depend on z' and are described by the term f z i.k 

/z',fc(A) « Ji(X) + 0{n 3 ), z e Ck- (28) 

Here the omitted higher-order terms correspond to 
D( z, z') > 2. We note that the first-order term in 
(28) becomes greater than the entropic term Fk (25) for 
M -k> log 2 Qun 2 ). 

Based on the above properties of Jd(z,z') we look for 
the solution of Eq. (23) in the avoided-crossing region as 
a sum of two terms 

fhz = u z + u z , z £ Ck, 1 < k < M (29) 

Here Uk and u z are slowly- and fast- varying functions of 
z, respectively, that satisfy the following equations: 

u z = 2 ~ n T [$ (A - aVo)- 1 + 8( A)*] , (30) 

^ k H z — ^ ^ i^Z)(z,z , )^z / Hz' H” T£ 0 ^ ^ yJ]D(z, w)^w- 

z'^£o wG>Co 

We note that the amplitudes \h w here correspond to the 
ground state(s) |w) of the final Hamiltonian, H P |w) = 
£o|w) and the decomposition (29) applies only to ampli- 
tudes with z ^ jC 0 . The above system of equations is 
closed by adding the equation for given by (23) with 
z = w and (29) taken into account. 

Eq. (30) for u z is coupled to the rest of equations via 
relation between the cumulative terms $ = $ + $ + $ 0 - 
Here $ = ]C Z 8 z u z is a part of $ (23) contributed by 
a slow- varying component and similarly, $ = 

We compute the cumulative term $ using equations for 
u z in (30) and also Eqs. (23)- (26) 

$ = -k(tiis( A)) $ 0 , ( 32 ) 

we£o 


where n(x) = x/(l + x ). In the initial stage of the al- 
gorithm the amplitudes of the solution states are small 
and |$ 0 | ~ |$| ~ 2 _ ". When these terms are neglected 
$ ps $, and (30) gives a closed-form algebraic equation 
for A which solution Aq(t) « a(r)V o — 2 pr + 0(p?) ac- 
curately tracks the adiabatic ground state energy, from 
r = 0, up until the small vicinity of the avoided-crossing, 
r w t* (see below) where |$o| ~ 1- In this region (and 
also for r > r*) one needs to retain the terms oc $ 0 in 
Eq.(30) that couple it to the equations for u z and \I/ W . 

We express the components u z in terms of the ampli- 
tudes by inverting the equation for u z in (30) with 
sequential iterations. The result is substituted into the 
sum J2 z ^c 0 i/d( z,w)U z that enters equation for U/ w (cf. 
(23)). A key observation is that the dominant contribu- 
tion to this sum comes from the states z that are one 
bit-flip away from the states w: r = d( z, w) = 1 (cf. 
Eq.(28)). This can be seen from the fact that amplitudes 
u z corresponding to not very large r <C n, are scaled with 
n as and also from the form of J r above. The 

contribution to the sum from the rest of terms (r > 1) is 
0(n~ 2 ). Finally, after excluding u z from the equations 
for in (23) and using (30), (32) we obtain a closed- 
form equation for A. We give it below in the region of 
interest \t — 1/2 1 <C 1 

(A - Aq(t)) (A - A f 0 (r) = — n 2 2“"A 2 /4 (33) 

A w d^ 2 (l + pr In 2 + 0(p 2 )) . 

The form of Aq(t) is given above and Aq(t) ps ts o + 
1/2 is another branch of the eigenvalues that eventually 
approaches £o for r — > 1. 

The minimum gap between the two roots of (33) is 
9m\n — n A 2 -ra / 2 and it is achieved at the point r = r* 
where the branches Aq(t), Aq(t) cross each other, t* m 
1/2 + \og 2 (do/p)/(4n). Based on the above analysis one 
can estimate the matrix element ~ n in (14). Then 
the complexity of the quantum adiabatic algorithm is, 
do#oi/Smin ~ (ndo) _1 2". To logarithmic accuracy the 
results are in agreement with the numerical simulations 
results given above. 

It follows from the discussion above that eigenvalue 
branch Aq(t) corresponds to an extended state J2z^ z \ z ) 
with exponentially small individual amplitudes, \u z \ ~ 
2 ~n/ 2 (cf. (30)). On the other hand Ag(r) corresponds to 
a localized state (|$o| ; |^| ~ 1) and therefore the coupling 
between the two states is exponentially small (cf. r.h.s. 
in Eqs. (30) and (33)). A similar picture applies to the 
avoided- crossing of the extended state energy Aq (r) with 
energies of localized states Aji. (r) corresponding to z e £* 
with k <C n. The existence of the extended eigenstate 
of H(t ) at later times, r > r*, whos properties do not 
depend on a particular instance of SPP can be explained 
as follows. According to (18)-(23), matrix elements of the 
Green function associated with the driver Hamiltonian 
(z'|(A — ay) -1 |z) depend only on a Hamming distance 
r = d( z,z'); on the other hand it follows from the above 
that the joint distribution of energies p r (E, E') factorizes 
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FIG. 1: Logarithmic plot of C m in vs n for randomly gener- 
ated instances of SPP with 25-bit precision numbers. Vertical 
sets of points indicate results of different trials (~ 100 trials 
for each n, except n— 17). Median values of C m in are shown 
with rectangles. Linear fit to the logarithmic plot of median 
values for n between 11 and 17 is shown by the line and gives 
In C m in ~ 0.55n (C m ; n ~ 2°' 8n ) . Very close result is obtained 
for the linear fit if all data points are used instead of the me- 
dian values. Insert: plot of C(T) vs T for n=15, precision 
b=25 bits, do=22. Point 1 indicated with the arrow refers to 
the minimal value of complexity at T — T m j n = 22.67 where 
the total population of a ground level po(Tmm) = 0.15. Point 
2 refers to the value of T where po(T) — 0.7. 


for not too large energies (cf. [18]), and therefore the 
structure of the problem is being averaged out during 
the quantum evolution that begins from the symmetric 
state (12). 

We also study the complexity of the algorithm by nu- 
merical integration of the time-dependent Schrodinger 
equation with Hamiltonian H(t) and initial state IV’(O)) 
defined in Eqs. (11) ,(12) ,(15)-(17) . Here we relax the 
condition % <C (E) used above in the analytical treat- 
ment of the problem; in simulations the value of M is set 
automatically to be an integer closest to log 2 Y^= o a i 
(cf. (16)). We introduce a complexity metric for the 
algorithm, C(T) = (1 + T)do/po(T) where Po(t) = 


Swe £ 0 lV’w(i)| 1 2 3 4 5 6 - A typical plot of C(T) for an instance 
of the problem with n=15 numbers is shown in the insert 
of Fig. 1. At very small T the wavefunction is close to 
the symmetric initial state and the complexity is ~ 2". 
The extremely sharp decrease in C(T) with T is due 
to the buildup of the population po(T) in the ground 
level, £ z = £o, as quantum evolution approaches the adi- 
abatic limit. At certain T = T m i n the function C(T) 
goes through the minimum: for T > T m i n the decrease in 
the number of trials do/po(T) does not compensate any- 
more for the overall increase in the runtime T for each 
trial. For a given problem instance the “minimal” com- 
plexity C m i n — C(T m i n ) is obtained via one dimensional 
minimization over T . Plot of the complexity C m i n for 
different values of n in Fig. 1 appears to indicate the 
exponential scaling law, C m i n ~ 2 0 - 8 " for not too small 
values of n > 11. 

In conclusion, we have studied both analytically and 
numerically, the complexity of the the quantum adiabatic 
algorithm applied to large random instances of SPP. We 
obtained the asymptotic behavior of the minimal gap 
that scales exponentially with problem size n. Analyt- 
ical results for the algorithm’s complexity are in qual- 
itative agreement with the results of numerical simula- 
tions of the algorithm for small-to-moderate values of n. 
We also developed a general technique for the analyzes 
of quantum spin-glass problems based on the separation 
of fast- and slow- varying parts of the wavefunction and 
the related picture of avoided crossing between the en- 
ergies of extended and localized states (as the system 
parameters vary). This approach can be applied to dif- 
ferent examples of random COPs where one searches for 
the solution using continuous-time quantum algorithms. 
Among the other applications is the analysis of tunneling 
phenomenon in the low-temperature dynamics of random 
magnets. 
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